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Wt  coniidtr  a  ■iaplifiad  satallitt  trajectory  problea, 
corresponding  to  a  flat  earth  assuiaptlon,  first  treated  by 
Okhotslmskll  and  Eneev.  We  present  a  numerical  solution 
based  upon  the  functional  equation  technique  of  dynamic 
programming,  and  a  proof  of  the  fundamental  result  In  the 
analytic  solution. 

The  same  computational  approach  can  be  applied  to 
more  realistic  trajectory  problems. 
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AN  APPLICATION  OP  DYNAMIC  PRDQRAMMINO 
TO  THE  DETERMINATION  OP 
OPTIMAL  SATELLITE  TRAJECTORIES 

Rlohard  Bellman 
Stuart  Dreyfus 


1,  Introduction 

In  this  paper  we  wish  to  consider  the  problem  of  deter¬ 
mining  optimal  satellite  trajectories.  To  illustrate  the 
genera),  techniques  of  dynamic  programming,  we  shall  consider 
In  detail  a  simplified  satellite  trajectory  problem,  posed 
and  treated  by  D.  E.  Olchotslmskll  and  T.  M.  Eneev  In  a 
paper  originally  published  In  RussleLn  In  [l]^»  and  then 
published  In  English  translation  In 

Despite  superficial  evidence  of  the  oorz* ••ectnese  of 
their  theories  and  results.  It  would  appear  that  the  mathe¬ 
matical  argumentation  In  their  papers  Is  at  best  Incomplete. 
One  purpose  of  this  paper  is  to  deduce  one  of  thel * 
principal  results  rigorously  using  standard  variational 
arguments.  What  Is  rather  Interesting  Is  that  we  arrive  at 
a  variational  problem  of  unconventional  type  which  has  not 
been  treated  to  any  extent  In  the  literature. 

The  principal  part  of  our  paper  Is  devoted  to  the  compu¬ 
tational  solution  of  the  problems  by  means  of  the  functional 
equation  technique  of  dynamic  programming.  We  shall  present 

*  [*1]  .  Uspekhl  Flzlohesklkh  Nauk,  September,  1957- 

••  L2]  •  J.  British  Interplanetary  Soc . ,  Jan.-Peb.,  1958. 
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some  numerical  results.  The  importance  of  this  approach 
lies  in  the  fact  that  we  can  solve  in  the  same  manner 
variational  probleme  which  defy  precise  analysis. 

2.  The  Simplified  Problem 

We  wish  to  ascertain  the  thrust  control  policy  and 
fuel  consumption  regime  which  will  put  a  satellite  into 
orbit  at  a  specified  altitude  with  maximum  horizontal 
Goniponent  of  velocity. 

Essential  almpllf Icatlons  arise  from  the  neglect 
of  aerodynamic  forces  and  the  assumption  that  the  terreetrlal 
gravitational  field  la  plane-parallel. 

Determination  of  paths  of  minimum  fuel,  maximum 
altitude,  and  so  on,  can  be  treated  along  the  same  lines  as 
the  following  discussion. 

3.  Mathematical  Fonculation 

The  equations  of  motion  of  a  satellite  traveling  over  a 
flat  earth  in  a  Cartesian  coordinate  syetem  will  be  taJcen 
to  be 

(1)  ^  =  p  cos 

^  =  p  sin  ^  -  g 


rovji 
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Her« 

(2)  (a)  X  and  y  are,  as  usual,  the  horizontal  and 
vertical  coordinates, 

(b)  u  and  w  are  the  horizontal  aind  vertical 
projectiona  of  velocity, 

(c)  p  is  the  magnitude  of  acceleration  due  to 
reaction  force, 

(d)  ^  is  the  inclination  of  the  thrust  to  the 
horizontal . 


If  we  Introduce  the  quantity  V  as  the  velocity  avail¬ 
able  to  the  satellite  in  the  idealized  case  of  no 
gravitational  force,  we  obtain  the  relation 


(3) 


dv 

IF 


P« 


The  variable  V  will  be  a  monotone  function  of  the  quantity 


of  fuel.  Since 
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(M 


where  P  is  the  thrust,  and 


(5) 


p  c  dM 


Where  M  !•  the  weight  end  o  is  the  exit  Telooitx  of  the 

gases,  we  cam  solve  for  M  in  terms  of  the  "Ideal  avail¬ 
able  velocity,"  V,  obtaining  the  equation 


(6) 


M  = 


V/c 


where  is  the  weight  of  the  empty  reoket. 

The  equations  of  motion,  (3.1 )»  together  with  (6)  pre¬ 
ceding,  which  yields  mass  as  a  function  of  V,  and  (4) 
above,  giving  acceleration  In  terms  of  thrust  and  mass, 
enable  us  to  determine  optimal  Inollnatlons  and  optimal 
magnitude  of  thrust  as  functions  of  V.  In  the  next  section 
we  shall  consider  the  associated  variational  problem. 

4.  Variational  Formulation 

A  dlrvfct  statement  of  the  problem  of  determining 
optimal  thrust  and  Inclination  leads  to  Immediate  diffi¬ 
culties  because  of  the  linearity  of  the  equations  and 
criterion  function.  An  optimal  policy  would  consist  of 
either  zero  or  Infinite  accelerations.  Since  this  last  Is 
physically  meaningless,  in  order  to  pose  a  sensible 
variational  problem  we  Impose  a  constraint  of  maximum 
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posiibla  thrust,  Q.  Purthermore,  in  these  initial  sections 
we  shall  agree  to  bum  all  the  fuel  at  nuucinoun  allowable 
rate  and  then  coast. 

Ttiis  assumption  will  not  be  made  in  later  sections 
where  we  obtain  a  direct  computational  solution  using 
dynamic  programming  techniques.  As  will  be  seen,  even  under 
these  simplifying  assumptions,  we  are  led  to  variational 
problems  of  some  novelty  and  interest. 

With  the  foregoing  assumptions,  the  problem  posed 
verbally  in  ^2  .ecoines  that  of  determining;  the  inclination 
function  fi{y)  which  maximizes 

(1)  J(^)  -  - oo«  ^  dV, 

^0 


subject  to  the  restriction  that 


and  the  relations  in  (3.1)  snd  (3.3)  -  (3.6).  Here 

(3)  (s)  H  is  the  prescribed  altitude, 

(b)  Vq  is  the  initial  ideal  velocity. 

The  functional  in  (2)  evaluates  the  altitude  gained 
during  burning  together  with  the  altitude  obtained  after 
burnout  due  to  the  vertical  component  of  velocity  at 


burnout. 


y 


Coasting 


Burnout  point 
Pig.  2 


3.  Preliminary  Tranafonnal.ion 

Referring  to  the  equations  of  (3>i)»  we  see  that  the 
constraint  in  (4.2)  may  be  written  in  the  form 


(1) 


H. 


Turning  back  to  (3*1)/  we  may  write 

(2)  w(V)  -  Wq  (.in  ^-£)dV. 

^0 


Using  this  relation  in  (5-i)»  we  obtain  the  relation 


where  the  function  r(V)  is  given  by 
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g  “  ir 


-  k,  +  k, 


V/c 


i^ert  k^  and  k^  az^  constants. 

6.  Diaouaaion 

V«  thus  obtain  the  problem  of  maximizing  J{^),  as 
glT«n  by  (4.1),  subject  to  the  constraint  of  (3)  above. 
This  la  rapreaentatlve  of  an  Interesting  class  of  varia¬ 
tional  problema  which  do  noc  appear  to  have  been  discussed 
in  any  detail  heretofore. 

TTie  general  problem  would  be  that  of  maximizing  a 
functional 

(1)  J^(y)  Pj^(x,y)dt, 


subject  to  constraints  of  the  form 


(2) 


JgCy)  -  P^(x,y)dt,  ..., 


P„(x,y)(ltj  -  k. 


We  shall  present  a  formal  analysis,  postponing  a  more 
rigorous  discussion  until  a  later  date. 


7«  Variational  Analysis 
As  us\ial«  we  set 
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(1)  + 


where  e  ie  an  infinitesimal.  We  have 


and  (5.3)  yields 

(3)  cos  7  r(V)dV  +  4  ?  cfvl  T COB  ?  dV  »  0. 

%  iW  i 


Since  is  arbitrary,  it  follows  that  there  exists 
a  constant  ^  such  that 


it  follows  that  we  obtain  the  important  conclusion  that  the 
optimal  policy,  is  characterized  by  the  property  that 

(6)  tan  ?  = 

a  constant . 

This  agrees  with  the  result  claimed  in  ^2j . 
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8,  Dyn&inic  Programming  Approach — I 

L«t  U8  now  see  how  we  oan  employ  the  functional  equation 
Approach  of  dynunic  prograosaing  to  obtain  a  computational 
solution.  The  basic  idea  is  to  regard  the  problema  of  the 
oaloulus  of  variations  as  particular  examples  of  multi-atage 
decision  processes  of  continuous  type.  This  approach  is 
discussed  in  some  detail  in  [3]  ••• ,  and  applied  to  optimal 
trajectory  problems  In 

The  state  variables  are  altitude  y,  vertical  component 
of  velocity  w,  and  ideal  available  velocity  V.  Conse¬ 
quently,  we  introduce  the  function 

(1)  f(V,w,y)  o  the  additional  horizontal  velocity  ob¬ 

tained  starting  at  altitude  y,  vertical 
component  of  velocity  w  and  ideal 
available  velocity  V,  and  using  an 
optimal  policy. 

Referring  to  the  defining  equations  of  motion  in  ^3,  and 
using  the  Principle  of  Optimality  as  in  [3]  or  pi)  ,  we  obtain 
the  functional  equation 

lij  .  Bellman,  R. ,  Pyemic  Programming,  Princeton 
university  Press,  Princeton,  Wew  Jersey, 

Ofl  •  Cartaino,  T.,  and  S.  Dreyfus,  "Application  of 
djnaalo  prograssolng  to  the  airplane  minimum  time- to-climb 
problen, "  Aeronautical  Engineering  Review,  vol.  16,  no.  6, 

1957,  PP.  tptt: -  - - 
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(2) 


f(V,w,y)  - 


Max 


Max 


cos  f(v  -  AV,w  4  ^V,y  4 


cos  ^AV  4  f  V  - 


«  e''/® 

iV,w  +  I  iln  fl - ip - 


y  +  -|jr—  iV) 


whers  AV  Is  regardsd  as  a  swill  quantity. 

Lstting  V  assume  only  a  finite  set  of  values 
0, AV,2AV, . . . ,NAV,  we  see  that  the  computation  becomes  that  of 
determining  a  sequence  of  functions  of  two  variables 
*  f(NA,w,y),  using  (2). 

9.  Dynasiic  Programming  Approach — II 

In  order  to  simplify  the  computation,  we  use  a  lAgrange 
multiplier  formalism,  as  discussed  in  [5] to  reduce  the 
problem  to  one  of  determining  a  sequence  of  functions  of  one 
variable . 

In  place  of  maximizing  J(^)  subject  to  the  constraint 
of  (5*1 )i  we  consider  the  problem  of  maximizing 

(1)  -y’®  ooB  ^(dv  ♦y’®  ^  dv, 

^0  ^0 


subject  to  the  constraints  of  the  equations  of  motion.  Here 
A  is  the  Lagrange  parameter. 


f5j  .  Bellman,  R.,  "Dynamic  programming  and 
Lagrange  multipliers,"  Proc.  Nat.  Acad.  Sol.  USA,  vol.  42, 
1956,  pp.  767-769. 
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The  new  fvmctional  equation  for  the  maximum  value  Is 


(2) 


f(V,w)  =  Max 


Max 


/wM  • 

COB  -f  -  * 


V/c 


W 


-AV 


M.e 


V/c 


+  f  V  -  AV,w  sin  ^ - - 


/wAt  +  f ( V,w  —  gAt ) 

where  the  second  alternative  within  the  Max  | 
represent!  a  decision  to  coast  for  a  small  time  Interval  At. 

The  parameter  A  la  adjusted  until  the  altitude  con¬ 
straint  of  (  4,2  )  Is  met.  By  using  the  Lagrange  parajneter, 

we  have  partitioned  a  pi’oblem  orlglnially  involving  a  sequence 
of  functions  of  two  variables  Into  a  set  of  problems  involving 
functions  of  one  variable.  The  gain  in  computing  time  emd 
effort  Is  considerable. 


10.  Computational  Aspects 

The  numerical  solution  is  obtained  by  iterating  the 
recurrence  equation  (  'J.2)  backwards  from  the  known  final 
values.  The  calculation  is  begun  by  observing  that,  if  bum 
out  occurs  with  a  vertical  component  of  velocity  w,  the 
additional  altitude  obtained  during  coasting  will  be  w  /2g 
and  the  additional  horizontal  velocity  will  be  zero.  Hence 
f(0,w)  -  Aw  /2g.  A  table  containing  f(0,w)  for  a  range  of 
w  values  (we  do  not  yet  know  to  what  burnout  value  of  w 
the  optimal  policy  will  lead)  Is  stored  in  the  high  speed 
memory  of  the  computer. 
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Thls  tabular  function  la  now  us«d  to  determint  a  new 
function,  f(AV,w),  the  total  additional  horizontal  velocity 
plus  A  times  the  altitude  that  can  be  attained  starting  with 
a  small  quantity  AV  of  "aYaJllable  velocity  (fuel)“  and 
vertical  velocity  component  \f.  This  calculation  Is  per¬ 
formed  using  equation  (4.2).  We  'actually  evaluate  the  gain 
associated  with  choices  of  different  f^'o  and  P's  and 
compare  this  with  the  return  from  a  decision  to  coast.  On 
this  basis  we  pick  the  optimal  decision.  The  return  from 
this  decision  Is  recorded  In  the  computer  memory  as  the  value 
of  f(AV,w)  for  the  particular  value  of  w  considered.  A 
second  table  Is  constructed  giving  the  optimal  decision  that 
yielded  f(AV,w).  A  third  table,  J(aV,w),  Is  maintained 
giving  the  total  altitude  gained  when  following  an  optimal 
path  starting  from  (aV,w).  Since  we  are  flying  so  as  to 
maximize  horizontal  veloclt.v  plus  A  times  altitude,  this 
third  table  Is  Just  a  convenient  record  that  Is  not  used  In 
the  calculation,  but  which,  when  the  Iteration  of  equation 
(4.2)  Is  finished,  yields  Immediately  the  total  altitude  (and 
hence  the  horizontal  velocity,  f(V,0)  —  AJ(V,0)),  gained  by 
following  an  optimal  trajectory. 

Once  the  technique  described  above  for  calculating 
f(AV,w)  using  the  table  of  f(0,w)  has  been  programmed  for 
a  computer,  it  is  a  simple  matter  to  have  the  same  code 
calculate  f(2AV,w)  from  f(AV,w)  and,  finally,  f(V,w) 
from  f(v  —  AV,w).  Notice  that  at  each  stageof  thin 
computation  only  one  table  of  the  function  f  Is  needed  to 
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oompute  the  next  table  in  the  sequence.  Once  a  table  has 
been  computed  and  used  in  the  calculation  of  the  next  table 
It  can  be  printed  by  tne  conyputer  and  destroyed  in  memory. 
Hence  the  computer  memory  capacity  required  is  determined  by 
the  number  of  discrete  points  chosen  for  the  w-table,  and 
does  not  depend  on  the  fineness  of  the  AV  grid.  The  total 
time  for  a  calculation  depends  inversely  on  the  sise  of  AV. 

At  the  completion  of  the  backwards  iteration  of  equation 
(4.2)  one  knows  the  horizontal  velocity  and  altitude  obtained 
by  an  optimal  policy  for  the  specified  initial  conditions. 
Also  ths  initial  decision  for  the  starting  point  is  deter¬ 
mined  by  the  nature  of  the  calculation  of  f(V,w).  To  re¬ 
construct  the  optimal  path  in  Its  entirety  one  now  determines 
the  new  value  of  w  after  using  the  prescribed  decision  for 
the  first  AV  interval.  In  calculating  f(V  —  AV,w)  for 
this  w- value,  an  optimal  decision  was  determined  and  lecoixied 
(since  the  actual  w  may  not  be  a  point  of  the  w-grld, inter¬ 
polation  may  be  necessairy)  and  this  decision  is  used  during 
the  Interval  V  —  AV  to  V  —  2a V,  In  this  manner  we  use  the 
output  of  the  sequence  of  calculations,  pixjcesslng  them  in 
the  opposite  order  frx)m  t^iat  in  which  they  wei^  computed. 

The  above  operation  may  be  performed  easily  by  the  com¬ 
puter  as  the  final  step  of  the  calculation  if  the  requisite 
tables  are  stored  on  tape  or  punched  into  cards. 

Once  the  problem  has  been  solved,  one  examines  the  final 
altitude  to  determine  if  the  required  height  was  attained.  A 
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new  value  of  the  Lagrange  multiplier  ^  la  then  calculated 
based  upon  pi'cvlous  value*  and  resul'  and  the  calculation 
la  repeated. 

One  calcul  yields  t>  >,tlmal  path,  In  terms  of 

horizontal  velocl'^  .  /A  -  aii.  iJe,  for  a  wide  range  of 
Initial  vertical  components  of  velocity.  This  variety  of 
results  Is  of  Interest  In  problems  where  Initial  vertical 
velocity  18  not  necessarily  specified  and  the  answers  for  a 
range  of  values  is  desired.  Secondly,  after  several  varia¬ 
tions  of  optimal  trajectories  to  several  different 

altitudes  are  Icnown,  yielding  an  Interesting  estimate  of  the 
trade-off  between  altitude  and  velocity  along  optimal 
trajectories , 

11.  Numerical  Results 

For  all  calculations,  we  have  assuii^ed  a  hypothetical 
missile  with  the  following  chara-terlstlcs : 

E^pty  weight,  -  3000  lbs. 

Exhaust  velocity,  c  -  11,000  ft ./sec. 

Maximum  thrust,  ^j^ax  "  30^^000  lbs. 

Minimum  thrust  with  engine  on,  ■  50*000  lbs. 

Total  Ideal  available  velocity  30,000  ft  ./sec. 

Tliese  data  ln^ly  a  total  weight  at  takeoff  of  76,43b 

lbs . 

A  value  of  A  of  .00142  yielded  a  final  altitude  of 
approximately  450  miles  with  a  horizontal  component  of 
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velocity  at  thlB  altitude  of  26,300  ft. /sec. 

Various  parameters  and  grid-sizes  required  for  numerical 
solution  were  chosen  as  follows: 

1.  AV  =  1000.  Therefore  the  recurrence  relation 
was  Iterated  30  times. 

2.  Aw  =  50.  Each  table  of  f(V,w)  contained 
281  numbers,  since  w  was  allowed  to  assume  value 
from  0  to  14,000. 

3*  A^  «  .01  mdlAn.  AdsilBslble  thrust  angles 

were  0,  .01,  .02,  ...,  t/2  radians. 

4.  Thrust  could  assume  values  300,000,  250,000, 
200,000,  150,000,  100,000,  or  50,000. 

These  numbers  were  determined  experimentally.  They 
possess  thf  property  that  a  further  refinement  has  little  or 
no  effect  on  the  computed  solution. 

A  condensed  summary  of  the  solution,  as  computed  on  the 
RAND  JOHNNIAC  computer  In  20  minutes.  Is  shown  below. 

It  should  be  noted  that,  although  In  this  simplified 
study  the  rocket  Is  flown  at  maximum  thirust  until  burnout  and 
the  thrust  direction  obeys  a  simple  j.aw,  the  computational 
scheme  assumes  neither  of  these  results.  It  Is  therefore 
applicable  to  more  general  problems  that  are  not  amenable  to 
conventional  mathematical  analysis. 
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V 

(ft/ 

sec  ) 

Hase 

(lb) 

h 

(ft) 

w 

(ft/ 

■*c ) 

u 

(ft/ 

oec ) 

(rad.  ) 

p 

(Iba) 

Time 

( sec ) 

30,000 

76,456 

0 

0 

0 

.560 

300,000 

0 

23 , 000 

48,529 

16,045 

1,514 

4,279 

.523 

300,000 

33.3 

20,000 

30,803 

37,720 

3,306 

8,625 

.501 

300,000 

54.4 

15 >000 

19.532 

10u,344 

5,239 

13,020 

.490 

300,000 

67.8 

10,000 

12,410 

152,954 

7,330 

17,436 

.480 

300,000 

76.3 

5,000 

7,877 

192,848 

9,464 

21,871 

.480 

300,000 

81.7 

burnout 

end  of 

5,000 

224,920 

11,650 

26,313 

.472 

300,000 

85.1 

coast 

5,000 

2,337,679 

0 

26,313 

0 

0 

444.7 

12.  Flow-chart 


A  flow  diagram  of  the  program  Is  shown  below: 
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.0  Y  J 


Yes 


Output 

Block  Transfer 
V  +  AV  - -  V 


No 
Yei 

_ ^ - 1  STOP 


y. 


